Attracting and repelling Lagrangian coherent structures 
from a single computation* 

Mohammad Farazmand 1 ' 2 and George Haller 2 
1 Department of Mathematics 
institute for Mechanical Systems 
ETH Zurich, 8092 Zurich, Switzerland 

January 22, 2013 



Abstract 

Hyperbolic Lagrangian Coherent Structures (LCSs) are locally most repelling or most at- 
tracting material surfaces in a finite-time dynamical system. To identify both types of hyperbolic 
LCSs at the same time instance, the standard practice has been to compute repelling LCSs from 
future data and attracting LCSs from past data. This approach tacitly assumes that coherent 
structures in the flow are fundamentally recurrent, and hence gives inconsistent results for tem- 
porally aperiodic systems. Here we resolve this inconsistency by showing how both repelling 
and attracting LCSs are computable at the same time instance from a single forward or a sin- 
gle backward run. These LCSs are obtained as surfaces normal to the weakest and strongest 
eigenvectors of the Cauchy-Green strain tensor. 

Repelling and attracting Lagrangian coherent structures (LCSs) are ma- 
terial surfaces that govern mixing patterns in complex dynamical systems. 
Recent developments made the accurate computation of both types of struc- 
tures possible, but not for the same data set: repelling LCSs are invariably 
obtained from future data, and attracting LCSs from past data. For tempo- 
rally aperiodic flows, this practice locates repelling and attracting LCSs for 
two different finite-time dynamical systems. Here we resolve this inconsis- 
tency by showing that both types of LCSs can be computed at the same time 
instance from the same data set. 

1 Introduction 

The differential equations governing a number of physical processes are only known as observa- 
tional or numerical data sets. Examples include oceanic and atmospheric particle motion, whose 
velocity field is only known at discrete locations, evolving aperiodically over a finite time-interval 
of availability. For such temporally aperiodic data sets, classic dynamical concepts-such as fixed 
points, periodic orbits, stable and unstable manifolds or chaotic attractors-are either undefined or 
nongeneric. 

Instead of relying on classic concepts, one may seek influential surfaces responsible for the forma- 
tion of observed trajectory patterns over a finite time frame of interest. Such a surface is necessarily 
a material surface, i.e., a codimension-one set of initial conditions evolving with the flow. Among 
material surfaces, an attracting Lagrangian Coherent Structure (LCS) is defined as a locally most 
attracting material surface in the phase space (Haller and Yuan, 2000; Haller, 2011). Repelling 
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LCSs are defined as locally most repelling material surfaces, i.e., attracting LCSs in backward-time. 
Repelling and attracting LCSs together are referred to as hyperbolic LCSs. Both heuristic detection 
methods (Peacock and Dabiri, 2010) and rigorous variational algorithms (Haller, 2011; Farazmand 
and Haller, 2012; Haller and Beron-Vera, 2012) are now available for their extraction from flow data. 

All available hyperbolic LCS methods fundamentally seek locations of large particle separation. 
They will highlight repelling LCS positions at some initial time t = a from a forward-time analysis 
of the flow over a finite time-interval [a, b]. Similarly, these methods reveal attracting LCSs at the 
final time t = b from a backward-time analysis of the flow over [a, b]. The complete hyperbolic LCS 
distribution at a fixed time t G [a, b] is, therefore, not directly available. 

Two main approaches have been employed to resolve this issue (see figure 1 for an illustration): 

1. Approach I: Divide the finite time interval of interest as [a, b] = [a, to] U [to , 6] . Compute 
repelling LCSs from a forward run over [to, 6], and attracting LCSs from the backward run 
over [a, to] (see, e.g., Lekien and Ross (2010); Lipinski and Mohseni (2010)). Both repelling 
and attracting LCSs are then obtained at the same time slice to- However, they correspond to 
two different finite-time dynamical systems: one defined over [a, to] and the other over [to, b]. 
This approach works well for a roughly T-periodic system, when to — a and b — to are integer 
multiples of T. In general, however, hyperbolic LCSs computed over [a, to] and over [to, b] do 
not evolve into each other as to is varied, and hence the resulting structures are not dynamically 
consistent. In addition, one cannot identify attracting LCSs at time a or repelling LCSs at 
time b from this approach. 

2. Approach II: Extract repelling LCSs at the initial time a from a forward run over [a, &]; extract 
attracting LCSs at the final time b from a backward run over [a, b]. Obtain repelling LCSs at 
any time to G [a, b] by advecting repelling LCSs from a to to under the flow. Similarly, obtain 
attracting LCSs at any time to G [a, b] by advecting attracting LCSs from b to to under the 
flow. This approach identifies LCSs based on the full available data, and provides dynamically 
consistent surfaces that evolve into each other as to varies (Haller, 2011; Farazmand and Haller, 
2012). Since the forward-time advection of a repelling LCS (as well as the backward-time 
advection of an attracting LCS) is numerically unstable (see figure 2), this approach requires 
extra care to suppress growing instabilities (Farazmand and Haller, 2012). Even under well- 
controlled instabilities, however, a further issue arises in near-incompressible flows: repelling 
LCSs shrink exponentially under forward-advection, and attracting LCSs shrink exponentially 
under backward-advection. Therefore, while the LCSs obtained in this fashion are dynamically 
consistent, they require substantial numerical effort to extract and may still reveal little about 
the dynamics. 

Here we develop a new approach that keeps the dynamical consistency of Approach II but 
eliminates the instability and shrinkage of advected LCSs. Our key observation is that attracting 
LCSs can also be recovered as codimension-one hypersurfaces normal to the weakest eigenvector 
field of the forward Cauchy-Green strain tensor. These stretch- surfaces are obtained from the same 
forward-time calculation that reveals repelling LCSs as strain- surfaces, i.e., codimension-one surfaces 
normal to the dominant eigenvector of the forward Cauchy-Green strain tensor (Farazmand and 
Haller, 2012). The locally most compressing strain-surfaces and the locally most expanding stretch- 
surfaces then reveal repelling and attracting LCSs at the same initial time a based on a single 
forward-time calculation over [a, b). 

We demonstrate the results on three examples: an autonomous Duffing oscillator (§5.1), a direct 
numerical simulation of two-dimensional turbulence (§5.2) and the three-dimensional classic ABC 
flow (§5.3). 
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Figure 1: Schematic illustration of Approach I (a) and Approach II (b) in the extended phase space. 

flow geometry 



repelling LCS 
at time a 




computed 
repelling LCS 



repelling LCS at time to 

B 



t = a 




C numerical advection 

1 of the computed LCS 



t = to > a 



Figure 2: The errors in the computation of a repelling LCS grow exponentially as the LCS is advected 
forwards in time. The same statement holds for the backward-time advection of an attracting LCS. 



2 Preliminaries and notation 

Consider the dynamical system 

x = u{x,t), xG/7cIR n , tel=[a,b], (1) 
where u : U x / — » R n is a sufficiently smooth velocity field. For t$,t G /, define the flow map 



x i->- x(t;t ,x ), 



(2) 



as the unique one-to-one map that takes the initial condition xq to its time-t position x(t]to,xo) 
under system (1). 
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The forward Cauchy-Green strain tensor over the time interval / is defined in terms of the flow 
gradient Vi 7 ^ as 

Cf =(VF b ) T VF a h . (3) 

At each initial condition xq G U, the tensor C^(xo) is represented by a symmetric, positive definite, 
n x n matrix with an orthonormal set of eigenvectors {£jJ(#o)}i</c<n, an d with a corresponding set 
of eigenvalues {\{(xo)}i<k<n satisfying 

Cf(x )Z f k (x ) = \{(xo)Z f k (x ), k G {1, 2, ••• , n}, (4a) 

< \{(x ) < X f 2 (x ) < • • • < A{(x ). (4b) 

These invariants of the Cauchy-Green strain tensor characterize the deformation experienced by 
trajectories starting close to xq. If a unit sphere is placed at #o, its image under the linearized flow 
map Vi 7 ^ will be an ellipsoid whose principal axes align with the eigenvectors {£{(£o)}i<fc<n and 
have corresponding lengths {A{(xo)}i<fc<n- 

Similarly, the backward Cauchy-Green strain tensor over the time interval I is defined as 

C b = ( V F 6 a ) T VF 6 a . (5) 

Its eigenvalues {A^(xo)}i</c< n and orthonormal eigenvectors {££(#o)}i<fc<n satisfy similar properties 
as those in equation (4). Their geometric meaning is similar to that of the invariants of C*^, but in 
backward time. 



3 Repelling and attracting LCSs 

A repelling LCS over the time interval I is a codimension-one material surface that is pointwise more 
repelling over / than any nearby material surface. If 1Z(t) represents the time-t position of such an 
LCS, then the initial LCS position 1Z(a) must be everywhere orthogonal to the most-stretching 
eigenvector ^ of the forward Cauchy-Green strain tensor (Haller, 2011; Haller and Beron-Vera, 
2012). Specifically, we must have 

T Xa ll(a) ±&(x a ), (6) 

for any point x a G 1Z(a), where T Xa 1Z(a) denotes the tangent space of 11(a) at point x a . 

Similarly, an attracting LCS over the time interval / is a codimension-one material surface that 
is pointwise more attracting over / than any nearby material surface. If A(t) is the time-t position 
of an attracting LCS, its final position A(b) satisfies 

T Xb A(b) J_ e n (x b ), (7) 

for all points x^ G A(b). That is, the time-6 position of attracting LCS is everywhere orthogonal to 
the eigenvector £^ of the backward Cauchy-Green strain tensor C b . 

The relation (6) enables the construction of repelling LCS candidates at time t = a, while (7) 
enables the construction of attracting LCS candidates at the final time t = b (see, e.g., Farazmand 
and Haller (2012); Hadjighasem, Farazmand, and Haller (2012)). Since LCSs are constructed as 
material surfaces, they move with the flow. Therefore, LCS positions at an intermediate time 
to G [a, b] are, in principle, uniquely determined by their end-positions: 

K(to) = Fp{K(a)), A(to) = Ft°(A(b)). (8) 

As discussed in the introduction, however, using the advection formulae (8) leads to numerical 
instabilities. This is because the material surfaces involved are unstable in the time direction they 
are advected in. This instability can only be controlled by employing a high-end numerical integra- 
tor which refines the advected surface when large stretching develops. Even under high-precision 
advection, however, the end-result is an exponentially shrinking surface which only captures subsets 
of the most influential material surfaces. 
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4 Main result 



Here we present a direct method to identify both attracting and repelling LCSs at the same time 
instance, using the same finite time-interval. These surfaces, therefore, are based on the assessment 
of the same finite-time dynamical system, avoiding the dynamical inconsistency we reviewed for 
Approach I in the Introduction. 

In particular, we show that the initial position of an attracting LCS, A(a), is everywhere orthog- 
onal to the weakest eigenvector ^( of the tensor '. This, together with the orthogonality of the 
initial repelling LCS position 71(a) to the dominant eigenvector £^ of C*^, allows for the simultaneous 
construction of attracting and repelling LCSs at time t = a, utilizing the same time interval [a, b]. 
All this renders the computation of the backward Cauchy-Green strain tensor C b unnecessary. 

Definition 1 (Strain-surface). Let M(t) be an (n — l)-dimensional smooth material surface in L 7 ", 
evolving under the flow map over the time interval / = [a, b] as M(t) = Fl(M(a)). Denote the 
tangent space of M. at a point x G M. by T X A4. 

(i) A4(t) is called a forward strain- surface if A4(a) is everywhere normal to the eigenvector field 

T Xa M(a)±£*(x a ), Vx a eM(a). 

(ii) A4(t) is called a backward strain- surface if A4(b) is everywhere normal to the eigenvector field 
£ n , i.e., 

T Xb M(b)±£(x b ), Vx b eM(b). 

Strain-surfaces are generalizations of the strainlines introduced in Farazmand and Haller (2012) 
and Haller and Beron-Vera (2012) in the theory of hyperbolic LCSs for two-dimensional flows. By 
contrast, the stretch-surfaces appearing in the following definition have not yet been used even in 
two-dimensional LCS detection. 

Definition 2 (Stretch-surface). Let M(i) be an (n— l)-dimensional material surface as in definition 
1. 

(i) M(t) is called a forward stretch- surface if M(a) is everywhere normal to the eigenvector field 
£i , i.e., 

T Xa M(a) ±t((x a ), \/x a e M(a). 

(ii) M(t) is called a backward stretch- surface if M(b) is everywhere normal to the eigenvector field 
£, b \ i- e -? 

T X6 M(b)±£(x b ), Vx b €M(b). 

By definition, the local orientation of a forward strain-surface is known at the initial time t = a. 
The following theorem determines the local orientation of the same strain-surface at the final time 
t = 6, rendering the forward-advection of the surface unnecessary. The same theorem provides the 
local orientation of backward strain-surfaces at the initial time t = a (see figure 3 for an illustration). 

Theorem 1. 

(i) Forward strain-surfaces coincide with backward stretch-surfaces, 
(ii) Backward strain-surfaces coincide with forward stretch-surfaces 

Proof. See Appendix A. □ 

The following corollary summarizes the implications of Theorem 1, along with known results 
from Haller (2011) and Farazmand and Haller (2012). 
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(a) 




t = a t = b 



Figure 3: (a) A forward strain-surface evolves into a backward stretch-surface, (b) A forward 
stretch-surface evolves into a backward strain-surface. 



Corollary 1. Let lZ(t) and A(t) be, respectively, repelling and attracting LCSs of the dynamical 
system (1). Then the following hold: 

(i) A repelling LCS, 7Z(t), is a forward strain-surface, i.e., 1Z(a) is everywhere orthogonal to 
the eigenvector field £f. Furthermore, 7Z(t) is also a backward stretch- surface, i.e., 71(b) is 
everywhere orthogonal to the eigenvector field £f. 

(ii) An attracting LCS, A(t), is a forward stretch-surface, i.e., A(a) is everywhere orthogonal to 
the eigenvector field Furthermore, A(t) is also a backward strain-surface, i.e., A(b) is 
everywhere orthogonal to the eigenvector field 

Among other things, the above corollary enables the visualization of attracting and repelling 
LCSs simultaneously at the initial time t = a of a finite time-interval [a, b] over which the underlying 
dynamical system is known (see section §5 below for examples). This only requires the computa- 
tion of the forward-time Cauchy-Green strain tensor , rendering backward-time computations 
unnecessary. 

5 Examples 

Here we demonstrate the application of corollary 1 on three examples: the classic Duffing oscillator, 
a two-dimensional turbulence simulation, and the classic ABC flow. In the two-dimensional case 
(i.e., n = 2), we refer to strain- and stretch-surfaces as strainlines and stretchlines, respectively. 
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5.1 Duffing oscillator 



Here we show that even for a two-dimensional autonomous system, stretchlines and strainlines act 
as de facto stable and unstable manifolds over finite time intervals. Indeed, over such intervals, sets 
of initial conditions will be seen to follow stretchlines in forward time. Only asymptotically do these 
initial conditions align with the well-known classic unstable manifolds. 




Figure 4: Trajectories of system (9). The homoclinic orbits are shown in red. 
Consider the unforced and undamped Duffing oscillator 

xi x 2 , 

x 2 =^x 1 -x\, (9) 

whose Hamiltonian H(x\,X2) = \x\ — Ax\ + x\ is conserved along the trajectories (see figure 4). 
The hyperbolic fixed point (0,0) of the system admits two homoclinic orbits (shown in red), which 
coincide with the stable and unstable manifolds of the fixed point. 

By Definition 1, forward strainlines over a finite time interval are everywhere orthogonal to the 
eigenvector field £| °f the forward strain tensor . As a result, strainlines are trajectories of the 
autonomous ordinary differential equation (ODE) 

r , (s) = £{(r(s)), r(0)=r , (10) 

where r : s \-> r(s) denotes parametrization by arc- length. Similarly, forward stretchlines are 
trajectories of the ODE 

p'(*)=^(p(*)), p(0)=po, (11) 

with p : s 4 p(s) denoting an arclength-parametrization. Since we are interested in the de facto 
finite-time stable and unstable manifolds passing through the hyperbolic fixed point (0,0), we set 
r =Po = (0,0). 

We observe that as the integration time T increases, the unique strainline and the unique stretch- 
line through the origin converge to their asymptotic limits. Figure 5 shows the convergence of these 
curves around the hyperbolic fixed point (0,0). For integration times T > 2, the computed strain- 
lines and stretchlines are virtually indistinguishable from their asymptotic limits. Therefore, in the 
following, we fix the integration time T = b — a = 2 with a = and b = 2. 

Note that while the strainline is indistinguishable from the stable manifold, the stretchline differs 
from the unstable manifold (see figure 5c). Stretchlines as de facto finite-time unstable manifolds 
define the directions along which passive tracers are observed to stretch. To demonstrate this, in 
figure 6, three disks with radii 10 -3 , 5 x 10 -3 and 10 -2 are initially centered at the origin. For short 



7 



-0.08 -0.06 -0.04 -0.02 0.02 0.04 0.06 0.08 

X\ 

(a) 



CM 



0.08 




0.08- 


0.06 




0.06- 


0.04 




0.04- 


0.02 




0.02 







CI 

H 


-0.02 




-0.02a 


-0.04 




-0.04 - 


-0.06 




-0.06 - 


-0.08 




-0.08 - 



-0.04-0.02 0.02 0.04 

X\ 

(b) 




-0.04-0.02 0.02 0.04 

X\ 

(c) 



Figure 5: (a) Forward stretchline through the origin for three integration times T = 0.5 (— x — ), 
T = 1 (— □— ) and T = 2 (— o — ). (b) Forward strainline for the same integration times, as in panel 
(a), (c) The asymptotic position of the strainline (— o — ) and the stretchline (— o — ) compared to 
the classic stable and unstable manifolds (black). 



advection times, the tracers elongate in the direction of the stretchline, not the unstable manifold. 
Unlike the classic unstable manifold, stretchlines evolve in time and only become invariant when 
viewed in the extended phase space of the (x, t) variables. For longer advection times (not presented 
here), the stretchline converges to the unstable manifold and becomes virtually indistinguishable 
from it. 
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Figure 6: (a) Classical stable and unstable manifolds (black) are shown together with the stretchline 
through the origin (magenta). Three blobs of tracers with radii 10 -3 (blue), 5 x 10 -3 (yellow) and 
10 -2 (red) are centered at the origin. The tracers and the manifolds are then advected to time 
t = 0.1 (b) t = 0.2 (c) and t = 0.4 (d). Over the time interval [0,2], the stretchline is the de facto 
unstable manifold for spreading tracers. For larger advection times, this de facto unstable manifold 
practically coincides the classic unstable manifold of the origin 

5.2 Two-dimensional turbulence 

We consider a two-dimensional velocity field u : U x R + — » R 2 , obtained as a numerical solution of 
the Navier-Stokes equations 



d+u + u • Vu 



-Vp + vAu + /, 
V-iz = 0, 
u(x, 0) = uq(x). 



(12) 
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The domain U = [0, 2ir] x [0, 2ir] is periodic in both spatial directions. The non-dimensional viscosity 
v is equal to 10 -5 . The forcing / is random in phase and active over the wave numbers 3.5 < k < 4.5. 
The initial condition uo is the instantaneous velocity field of a decaying turbulent flow. We solve 
equations (12) by a standard pseudo-spectral method with 512 x 512 modes. The time integration 
is carried out by a 4th order Runge-Kutta method with adaptive step-size (MATLAB's ODE45). 
Equation (12) is solved over the time interval / = [0,50]. 

One can, in principle, compute an attracting LCS at the beginning of a time interval / = [a, b] 
by advecting the attracting LCS extracted at t = b back to t = a. As mentioned in the Introduction, 
however, this process is numerically unstable since attracting LCSs become unstable in backward 
time. Their instability is apparent in figure 7, where an attracting LCS (red) is advected backwards 
from t = 50 to the initial time t = 0. The advected curve is noisy and deviates from the true 
pre-image (blue curve). The true pre-image, the stretchline, is computed as a trajectory of the 
eigenvector filed £2 °f the forward Cauchy-Green strain tensor 




Figure 7: Stretchline (blue) and the advected image of an attracting LCS (red) at t = 0. The 
exponential growth of errors in backward-time advection of the LCS results in a jagged curve that 
deviates from the true attracting LCS. 

We now extract the set of attracting LCSs that shape observed global tracer patterns in this 
turbulent flow. Corollary 1 establishes that such LCSs are necessarily forward stretchlines, i..e, 
trajectories of (11). It then remains to select the trajectories of this ODE that stretch more under 
forward advection than any neighboring stretchline (Haller and Beron-Vera, 2012). 

The relative stretching of a material line is defined as the ratio of its length at the final time t = b 
to its initial length at time t = a. For a forward-time stretchline 7, one can show (see Appendix B) 
that the relative stretching is given by 




where £(7) is the length of 7 at time t = a. Note that no material line advection is required for 
computing the relative stretching in (13). 

In order to locate the stretchlines that locally maximize the relative stretching (13), we adopt the 
numerical procedure outlined in Haller and Beron-Vera (2012) for locating the locally least-stretching 
strainlines. Specifically, we first compute a dense enough set of stretchlines as the trajectories of 
ODE (11). We stop the integration once the stretchline reaches a singularity of the tensor field 
or crosses an elliptic transport barrier. 

A singularity of is a point where equals the identity tensor, and hence its eigenvectors 
are not uniquely defined (see Delmarcelle and Hesselink (1994) and Tricoche, Scheuermann, and 
Hagen (2000) for more details). An elliptic barrier is the outermost member of a nested set of closed 
curves that preserve their initial length (at time t = a) under advection up to time t = b (Haller 
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3. 




(d) 

Figure 8: (a) The concentric tracers with radii 0.05 (blue), 0.1 (yellow) and 0.2 (red). The stretchline 
(black) passing through the center is computed from the time interval [0, 50] (i.e., a = and b = 50). 
The tracers and the stretchline are then advected forward in time to t = 10 (b), t = 15 (c), t = 25 

(d). 

and Beron-Vera, 2012). In an incompressible flow, an elliptic barrier also preserves its enclosed area 
under advection, and hence the elliptic domain it encloses remains highly coherent. For this reason, 
elliptic barriers can be considered as generalizations of outermost KAM curves generically observed 
in temporally periodic two-dimensional flows (Haller and Beron-Vera, 2012). 

We locate elliptic barriers using the detection algorithm developed in Haller and Beron-Vera 
(2012) and Hadjighasem, Farazmand, and Haller (2012). With the location of these barriers and 
of the singularities of at hand, stretchlines are truncated to compact line segments, rendering 
the integral in (13) well-defined. Attracting LCSs at t = a are then located as stretchline segments 
that have higher relative stretching (13) than any of their C 1 -close neighbors. This process is briefly 
summarized in the following algorithm: 

Algorithm 1. 

1. Compute the Cauchy-Green strain tensor over a uniform grid. 

2. Locate elliptic barriers by the procedure described in Haller and Beron-Vera (2012) and 
Hadjighasem, Farazmand, and Haller (2012). 

3. Compute stretchlines as trajectories of (11). The initial conditions po are chosen from a uniform 
grid over the phase space. 
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Figure 9: (a) Forward stretchlines at t = 0. The attracting LCSs (i.e., locally most-stretching 
stretchlines) are highlighted in red. The green closed curves show the boundaries of elliptic regions. 
Tracers (blue circles) are used to visualize the overall mixing patterns, (b) Advected image of the 
attracting LCSs, tracers and elliptic barriers at time t = 50. 

4. Stop the stretchline integration once the stretchlines reach either a singular point or an elliptic 
region bounded by an elliptic barrier. 

5. For each stretchline so obtained, compute the relative stretching (13). 

6. Locate attracting LCSs as the stretchlines with locally maximal relative stretching. 

To illustrate the defining role of stretchlines in the formation of turbulent mixing patterns, we 
consider three concentric circles of tracers with radii 0.05, 0.1 and 0.2 at the initial time t = a = 
(see figure 8). The circles are centered on a stretchline with locally largest relative stretching (black 
curve). Then the stretchlines and tracers are advected to times to = 10, to = 15 and to = 25. In each 
case, we find that the tracer pattern stretches and alines with the evolving stretchline, as expected. 

We now turn to the global geometry of the attracting LCSs. Figure 9a shows stretchlines com- 
puted from a uniform grid of 30 x 30 points. Attracting LCSs at time t = 0, extracted as stretchlines 
with the locally largest relative stretching, are highlighted in red. Also shown are the elliptic barriers 
(greed closed curves), as well as a select set of blue tracer disks that will be used to illustrate the 
role of attracting LCSs. The advected positions of attracting LCSs, elliptic barriers and tracer disks 
are shown in figure 9b. Note how the attracting LCSs govern the deformation of the tracer disks 
in the turbulent mixing region. Meanwhile, the elliptic barriers keep their coherence by preserving 
their arclength and enclosed area. 

5.3 ABC flow 

In two dimensions, stretchlines are constructed as trajectories of the eigenvector field The 
resulting curves are, by construction, everywhere orthogonal to the eigenvector field £{. In higher 
dimensions, however, constructing stretch-surfaces that are everywhere orthogonal to the eigenvector 
is nontrivial. In fact, for a given eigenvector field, such a surface may only exists locally if a 
Frobenius condition is satisfied (Lee, 2009). This condition requires the eigenvectors spanning the 
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Figure 10: (a) A spherical tracer surface (blue) at time t = and the corresponding approximate 
stretch-surface (red) passing through its origin, (b) The advected positions of these surfaces at the 
final time t = 4. 

tangent space of the manifold (here, }2</c<n) to be in involution, i.e., their Poisson brackets 
[£/>£/] should be in the tangent space of the manifold for any z, j £ {2, 3, • • • , n}. 

Even when the subset of the phase space satisfying this Frobenius condition is known, construct- 
ing stretch-surfaces globally as smooth parametrized manifolds normal to a specific vector field is 
challenging (Palmerius, Cooper, and Ynnerman, 2009; Balzer, 2012). Here we only illustrate that 
locally constructed stretch-surfaces do govern the formation of tracer patterns in three-dimensional 
flows as well. 

We use the classic ABC flow (Arnold and Khesin, 1998) 

xi = Asm(xs) + C cos(x2), 
±2 = B sin(xi) + Acos(x3), 

xs = Csin(x2) + B cos(xi), (14) 

with A = 1, B = y^2/3 and C = y/T/3. The C f strain tensor is computed over the time interval 
/ = [0,4] (i.e., a = and 6 = 4). We release a spherical blob of initial conditions centered at (jr, tt) 
with radius 0.1. We approximate the stretch-surface passing through this point by the plane normal 
to the first eigenvector of . Figure 10a shows this plane together with the sphere of tracers 
at time t = 0. The advected images of the tracer and the plane at time t = 4 are shown in figure 
10b, demonstrating that the stretch-surface through the center of the tracer blob acts as a de facto 
unstable manifold in this three-dimensional example as well. 

6 Conclusions 

We have shown that both repelling and attracting LCSs (finite-time stable and unstable manifolds) 
at a time instance t = a can be extracted from a single forward-time computation over a time interval 
/ = [a, b]. This extraction requires the computation of the eigenvectors of the forward Cauchy-Green 
strain tensor C f . It has been found previously (Haller, 2011; Haller and Beron-Vera, 2012) that at 
time t = a, the position of repelling LCSs are strain-surfaces, i.e., are everywhere orthogonal to 
the dominant eigenvector of C$ . Here we proved that the t = a positions of attracting LCSs are 
stretch-surfaces, i.e., are everywhere orthogonal to the weakest eigenvector of . 

The attracting LCSs obtained in this fashion are observed as centerpieces around which tracer 
patterns develop. Even in autonomous dynamical systems, these evolving centerpieces of trajectory 
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evolution differ from classic unstable manifolds, forming de facto unstable manifolds over finite times. 

In two-dimensional dynamical systems, stretchlines can be directly computed as most-stretching 
trajectories of the autonomous ODE (11). In higher dimensions, stretch-surfaces satisfy linear 
systems of partial differential equations (PDEs), as any surface normal to a given vector field does 
(Palmerius, Cooper, and Ynnerman, 2009). While a self-consistent global solution of these PDEs 
remains numerically challenging, here we have illustrated the local organizing role of stretch-surfaces 
through the advection of their tangent spaces in the classic ABC flow. Results on the construction 
of attracting LCSs from globally computed stretch-surfaces will be reported elsewhere. 

Acknowledgements. G. H. acknowledges partial support by the Canadian NSERC under grant 
401839-11. 

Appendix A Proof of Theorem 1 

In order to prove Theorem 1, we need two lemmas. The first lemma draws a connection between 
eigenvalues of the forward- and backward-time Cauchy-Green strain tensors. The second lemma 
establishes a relation between their eigenvectors. 

Lemma 1. The largest eigenvalue of the forward-time strain tensor at a point x a G U 
coincides with the reciprocal of the smallest eigenvalue \\ of the backward-time strain tensor C b at 
the point = F^(x a ), i.e., 

Similarly, we have 

Proof. This follows directly from equation (13) in Haller and Sapsis (2011). □ 
Lemma 2. For any x a G Z7, the following identities hold for any k G {1, 2, • • • , n}. 

(&(x a ), VF?(x b )&(x b )) = \ f n (x a )\ b k (x b )(&(x a ), VF b a (x b )e k (x b )), (17) 

(e n (x b ),VF b a (x a )tl(x a )) = X b n (x b )X{(x a )(e n (x b )yF b a (x a )^ k (x a )), (18) 
where (•, •) is the Euclidean inner product between two vectors. 

Proof We prove identity (17). The proof of (18) is similar and will be omitted. 

First, note that since the flow map is invertible, we have (i^(x a )) = x a for any x a G U. 
Differentiating this identity with respect to x a , we obtain 

VF b a (x b )=[VF b (x a )]- 1 . (19) 
The result then follows from the identity 

(&M,VF b a (x b )e k (x b )) = (£f(x a ), [VF b a (x b )]- T [VF b a (x b )] T VF b a (x b )e k (x b )) 
= ([VF b a (x b )]- 1 ti(x a ),C b (x b )e k (x b )) 
= X b (x b )(\7F b (x a )e n (x a ),e k (x b )) 

= A^^)([VF a t ( a;o )]- T [VF a 6 ( a;o )] T VF 6 (a ;o )4 / ^),^fe)) 
= \ b k (x b )(Cf(x a )&(x a ), [VJjCca)]- 1 ^)) 
= \l(x a )\ b (x b )(&(x a ),VF b a (x b )e k (x b )), 

where we have used identity (19) twice. □ 
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Now we turn to the proof of Theorem 1. 
Proof of Theorem 1: 



(i) Assume that A4(t) is a backward stretch-surface. Then, by definition, A4(b) is everywhere 
orthogonal to the eigenvector field In order to show that A4(t) is a forward strain-surface, 
it suffices to show that M.(a) = Fg(A4(b)) is everywhere normal to the eigenvector field 
Since T Xh M(b) = span{££(#b)}2<fc<n for any £5 G M,(b), we have 

T Xa M(a) = span{VF 6 a (x 6 )^(x 6 )} 2 

for all x a := Fg(xb) G jVf(a). Therefore, it suffices to show that £l(x a ) _L VF^x^^^x^) for 
any x a G .M(a) and fc G {2, 3, • • • , n}. 

From Lemma 2, we have 

(£[(*„), VF b a (x b )^(x b )) = \l(x a )\Ux b )(ti(x a ), VF?(x b )&(x b )), (20) 

for any x a G A4(a) and fc G {2, 3, • • • , n}. 
Using identity (15), we obtain 

(ti(x a ),VF b a (x b )e k (x b )) = >^{Zf( Xa ),VF?(x b )&{x b )). (21) 

Hence, if 

\\(x b )±\\{x h ), fc€{2,3,--- ,n}, (22) 

then we have 

(^(x a ),VF 6 a (x fc )^(x fc ))=0, (23) 

for any fc G {2,3, •• • , n}. But since \\ < \\ < • • • < A^, conditions (22) hold if and only if 
^i(xb) 7^ ^2(^5). This condition holds away from repeated eigenvalues of C b . 

In short, if gftx b ) _L T Xh M(b) for all x 6 G then &(x a ) _L T Xa M(a) for any x a G .M(a) 

which implies that .M(a) is a forward strain-surface. This concludes the sufficiency condition 
of Theorem l-(i). 

As for the necessity of the same condition, let A4(t) be a forward strain-surface, i.e. T Xa M.(a) = 
span{^(x a )}i</e< n _i for any x a G M(a). Therefore, the tangent space of its advected image 
M.(b) is given by 

T Xb M(b) = S ptm{VF b a (x a )tl(x a )} Kk<n-1- 

To show that M(t) is a backward stretch-surface, it suffices to show that _L V F^{x a )^ k {x a ) 

for any x& G M(b) and fcG{l,2,---,n — 1}. Similarly to equation (21), one can show that 

VF h a (x a )e k (x a )) = ^4(^(x 6 ), VF h a {x a )e k {x a )), (24) 

which implies that VF^(x a )^l(x a )) = for fc G {1, 2, • • • , n — 1} away from the degen- 

erate points where A^ = A^_ x . 

(ii) The proof is identical to that of part (i). 

□ 
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Appendix B Relative stretching of stretchlines 

Here, we derive formula (13) for the relative stretching of forward stretchlines. Let jt be a smooth 
material line. Denote its time-a and time-6 positions by j a and 75, respectively. Then, the relative 
stretching of the material line j t over the time interval I = [a, b] is defined as 

9 (7,):=g, (25) 

where t denotes the length of a curve. 

Let r : s \-> r(s) be the parametrization of j a by arc- length, i.e., let \r'(s)\ — 1 for all s G [0, i(j a )]- 
Since 75 = F^(7 a ), the mapping F^or : s \-> F^(r(s)) is a parametrization of the curve 75. Therefore, 
its length ^(75) is given by 

«= r 7o) iv^(K«)y(«)id S 

./o 

= J y / (r'( S ),C/(r( S ))r'( S ))d S . (26) 

Now, if the material line j t is a forward stretchline, we have r'(s) = &(r(s)) for all s G [0,^)]. 
Substituting this in equation (26), we obtain 



£( lb ) = f {la) y/\£(r(s))d8 := f J^ds. 



Therefore, by definition (25), the relative stretching of a forward-time stretchline j t is given by 
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